Quantitative analysis of myofiber type composition in human and mouse skeletal muscles

Summary Skeletal muscles are composed of different myofiber types characterized by the expression of myosin heavy chain isoforms, which can be affected by physical activity, aging, and pathological conditions. Here, we present a step-by-step high-throughput semi-automated approach for performing myofiber type quantification of entire human or mouse muscle tissue sections, including immunofluorescence staining, image acquisition, processing, and quantification. For complete details on the use and execution of this protocol, please refer to Abbassi-Daloii et al. (2022).1


Cryosectioning of muscle biopsies
Timing: 20-30 min per muscle biopsy The purpose of this step is to cryosection muscle biopsies for immunofluorescence staining.
Note: For users unfamiliar with cryosectioning, we recommend reading the article by Ross et al., 10 which provides a detailed protocol.
Here, we describe the procedure for collecting cryosections: 1. Clean all the equipment with 70% alcohol. 2. Adjust the chamber temperature and the object temperature to À20 C and À22 C, respectively. 3. Place all the materials required for cutting (e.g., tweezers, specimen holder) inside the cryostat chamber to equilibrate to the temperature. 4. Install the anti-roll glass and a blade inside the cryostat chamber. 5. Adjust the thickness (we recommend cryosections of 16 mm thick for human tissue and 10 mm thick for mouse, but 8-16 mm thick leads to acceptable results). 6. Transport the muscle biopsy in liquid nitrogen to the cryostat. 7. Leave the biopsy inside the cryostat chamber to equilibrate to the temperature for 15-30 min. 8. Put some Tissue-Tek (depending on the biopsy size) on the specimen holder and place the muscle biopsy in the Tissue-Tek (using tweezers). 9. Place the specimen holder in the cryostat block when the muscle biopsy is completely fixed. 10. Collect the cryosections onto the SuperFrost slide.
Note: The tissue quality should be assessed by checking the first section under a normal bright-field microscope. Samples with extensive freezing damage should be excluded. In addition, since the purpose of this protocol is to perform myofiber typing for which cross-sections of myofibers are required, therefore you should examine if there are more cross-sections of myofibers rather than longitudinal sections of myofibers ( Figure 2). If needed, the muscle biopsy should be detached from the specimen holder and placed again with the correct orientation to get cross-sections. 11. Per biopsy, collect three cryosections onto SuperFrost slides.
Note: These three cryosections are used as technical replicates. Therefore, when pasting them onto a slide, make sure there is enough space around each cryosection to draw a circle using immunopen (step 15).
12. Store the slides at À20 C prior to staining.
Note: The slides can be stored at À20 C for several years.

Timing: A day and a half per experimental batch
This section describes immunofluorescence staining using antibodies for three myosin heavy chain (MyHC) isoforms (MyHC1, MyHC2A, and MyHC2X) and laminin. Our protocol for immunofluorescence staining with the antibody mix of MyHC isoforms and laminin was described. 11 Note: For conducting immunofluorescence staining in multiple experimental batches samples across groups should be randomized in different batches using a single antibody master mix. The antibody master mix can be stored at 4 degrees for 1-2 weeks.
Here, we provide step-by-step directions for conducting this experiment: 13. If slides were stored at À20 C, air dry them for 30 min at room temperature (RT). This incubation is important to limit detachment of the tissue during the incubations. 14. Outline each cryosection with an immunopen about 2-3 mm from the tissue edge.
Note: Do not draw the line too close to the muscle cryosections as it will introduce an artifact in the image processing step. 15. Wash the cryosections in PBST in the staining box. 16. Blocking.
17. Wash the slides three consecutive times with a large volume of PBST, each time for 5 min.

Primary antibody incubation.
Incubate sections with a mixture of the following primary antibodies for 2 h at RT: a. Rabbit anti-laminin, b. Mouse anti-6H1 detecting MyHC2X.
Note: The volume of the primary antibodies depends on the section area and you should make sure that the entire cryosection is covered with the antibodies.
19. Wash the slides three consecutive times with an excessive volume of PBST, each time for 5 min. Note: Make sure that the sections don't dry out overnight. The slides can be carefully placed on wet tissue.
23. Wash the slides three consecutive times with an excessive volume of PBST for 5 min. 24. Wash the slides once with an excessive volume of PBS for 5 min.
25. Mounting. a. Cover the sections with ProLongä Gold antifade reagent. b. Cover the slide with a cover slip.
Note: Avoid any air bubbles on the sections as they will affect the image acquisition.
c. Fix the cover slip with nail polish.
26. Place the mounted sample on a flat, dry surface. 27. Incubate for 24 h at room temperature in the dark. 28. Store slides at 4 C prior to imaging.

Image acquisition
Timing: 10-15 min (per sample) Here, we describe the imaging of the entire muscle sections using an Axio Scan.Z1 slidescanner (Carl Zeiss, Germany) image capturing using ZEN 2 (blue edition) software (the v2.6 was used in this protocol).
Alternatives: Other microscopes with high capacity image acquisition equipped with four fluorescence channels can be used.
Note: The imaging settings should be optimized on a test slide to specify the exposure time and intensity per fluorophore, as exposure time and focusing algorithm fade the fluorophore signal.
Note: We recommend using the channel with the highest signal-to-noise ratio (in this dataset the MyHC2A conjugated-Alexa Fluorâ 594 channel) to define the focus plane.
Note: For each fluorophore/channel, the intensity and exposure time should be optimized to get the best signal-to-noise ratio, without bleaching the fluorophores.
29. Make images with a 103/0.45 Plan-Apochromat objective lens ( Figure 3). Note: The same image acquisition settings must be used for all slides over all batches.
Note: When you perform the staining in multiple batches, the imaging should be also done in the same order and batches to keep the same duration between staining and imaging for all batches.
Note: For each slide, the output is a Carl Zeiss Image format (CZI) dataset, which includes an image for each section.

Image preprocessing
Timing: 1-2 min (per sample) Note: A video walkthrough of the image preprocessing and image processing steps described below can be found in Methods video S1.
Before the acquired images can be processed, the data needs to be curated, and some preprocessing might be required as explained below. The result of steps 32 and 33 should be a folder of multichannel tiff files with a pixel size between 1 and 5 mm that can be imported by Fiji. 6 31. Calculating shading profile. If the acquired images exhibit a significant amount of shading ( Figure 4A), this should be corrected since these intensity variations are not related to differences between distinct myofibers ( Figure 3). The acquisition software offers a number of functions for shading correction: as an online process during the acquisition, or as a post-processing step.

OPEN ACCESS
When these options do not give satisfactory results, it is possible to calculate 'averaged' shading profiles and use these for post-processing. When imaging with Axio Scan.Z1 slidescanner, this is achieved using ZEN Lite (v3.3 used for this protocol) and Fiji as below: a. Calculate a shading profile using the 'Shading Reference From Tile Image' in ZEN Lite for each channel in each slide.
i. This will produce a shading profile for each channel per slide. b. Load all shading profiles for one specific channel into Fiji and combine them into a stack using the 'Images to stack' command. c. Using 'Z Project' with 'Projection type' set to 'Median', calculate the median shading profile for this channel. d. Repeat steps b-c for all channels. e. Use the median shading profiles for the shading correction using 'Shading Correction' in ZEN Lite. This can only be done one channel at a time ( Figure 4).

Converting image format.
Depending on the slide scanner used, this can be achieved in a multitude of ways. When using the Axio Scan.Z1 slidescanner the following procedure can be used in Fiji: a. Open Fiji. b. Run ''0.Convert_CZI_to_Tiff.ijm'' macro to convert the slidescanner datasets from CZI to multichannel 16-bit TIFF files using BioFormats. 12 c. Provide a directory with CZI datasets.
i. Running this macro, the images will be downsampled (43) by averaging to improve the processing speed and reduce the required data storage. ii. The effective pixel size is 2.6 mm after downsampling.
Note: The effective pixel size will depend on the magnification of the acquisition. One should modify this to achieve a pixel size between 1 and 5 mm. It can be tuned by changing the scale parameter in ''0.Convert_CZI_to_Tiff.ijm'' macro (line 39). The amount of downscaling should be evaluated carefully by checking the output of the subsequent segmentation steps. If the segmentation quality is not sufficient, consider reducing the amount of downsampling.
d. For each image in each CZI dataset, a tiff file that ends with ''_s[X]_merged'' will be saved in the input directory provided.
i. X shows the image number in the CZI dataset (starting from zero).

Image processing
Timing: 5-10 min (per sample) Image processing will be performed using Fiji and Ilastik 7 software.
Note: A modular set of macros that process each step independently is created. All the macros are publicly available on GitHub: https://github.com/tabbassidaloii/ImageProcessing/tree/ main/MyofiberTyping/Macros. In the steps below, we specify which macro should be used.
Note: Macros save intermediate image files to facilitate debugging and rerunning certain steps when required.
33. Generating tissue mask. This step is a semi-automated process that uses an automated mask generation algorithm, followed by a manual step to check and correct the generated masks.
Note: The aim of these manual corrections is to remove artifacts such as tissue folds, out-offocus regions, scratches, and dirt objects.
Note: The image processing steps below expect a folder that contains multichannel 16  i. The first image will pop up with a mask generated automatically and a question box: ''Is the mask OK?'' ii. Press the OK button, if the generated mask is accurate. Otherwise, the mask can be adjusted ( Figure 5): Zoom in or zoom out on the area of interest by holding Ctrl and scrolling the mouse wheel.

OPEN ACCESS
Hold Shift and Right Click on the area of the image that needed to be added to the mask; Hold Alt and Right Click on the area of the image that needed to be removed from the mask.
Press OK when the mask is accurate to go to the next image.
Note: This manual process is the most time-consuming step in image processing.
e. For each image, a tiff file that ends with ''_Mask'' will be saved in the input directory provided.
34. Generating 'masked' copy of the laminin channel. In this step, to reduce any possible artifacts due to this binary mask, a gaussian blur is applied to the mask and the pixel values of the laminin channel outside the mask are set to the median intensity of these pixels.  Select the pen with the proper size to annotate the pixels which belong to the class selected.
Select the other class and repeat step 3. Evaluate the classifier performance by clicking on ''Live Update'' button and selecting ''Segmentation'' option (from ''Group Visibility'').
If required, improve the pixel annotation in each class. Train and evaluate the classifier on at least three different input images by selecting them from ''Current View''.
Note: Annotate pixels in multiple areas of each section representing the tissue characteristics.
Note: Avoid overtraining the algorithm on a single section or area, as this will reduce the classifier performance across the whole dataset.
Note: When the classifier works as you expect, we recommend continuing with image processing and evaluating the performance of the classifier based on the output of the next step.
i. Select ''4. Prediction Export'' from the menu on the left.
ii. Define the export image settings by clicking on ''Choose Export Image Settings'' as shown in Figure 7. iii. Save the project.
Note: An example of an Ilastik classifier file can be found on GitHub (https://github.com/ tabbassidaloii/ImageProcessing/tree/main/MyofiberTyping/Macros) under the name of '3.Pixelclass_Lamin_Masked.ilp'. i. Use this classifier to process all images. Open Command Prompt. Run the command below by providing paths required: For each image, a tiff file that ends with ''Masked_Probabilities'' will be saved in the image directory.

Laminin segmentation and myofiber region-of-interest (ROI) generation.
In this step, laminin segmentation is used to generate the regions-of-interest (ROI) (individual myofibers) for each image. a. Open Fiji. b. Run ''4.Segment_Lamin.ijm'' macro to segment and generate the ROIs for each image. c. Provide a directory with ''Masked_Probabilities'' tiff images. d. For each image. i. a tiff file (with ''_Segmentation'' extension) in the input directory provided. and, ii. an ROI file that ends with ''_ROI'' in the ROI subdirectory will be saved ( Figure 8).

Extracting the mean-fluorescence-intensity and ROI properties.
In  is high when the pixel-classification is high for the 'myofiber boundary' class all around the myofiber and low in the interior of the myofiber.
a. Open Fiji. b. Run ''5.Export_MFI_and_Laminin_Int_and_Distance.ijm'' macro for the laminin segmentation and generate the regions-of-interest (ROI) (individual myofibers) for each image. c. Provide a directory with tiff images. d. For each image, a file and three images will be saved.
i. a tab-delimited text file that ends with ''_MFI'' in the ROI subdirectory (Table 2).
ii. three jpeg images that end with ''_check[X]'' in the check subdirectory.
Note: The distance of each myofiber to the edge of the tissue (Table 2) can be used to remove myofibers in case of suspected artifacts, such as close to the tissue borders (troubleshooting 2).
Note: To aid interpretation and verification of the results, visualizations of the segmentation, as well as the measured parameters such as mean fluorescence intensity are generated (jpeg images ending with ''_check[X]'').

Myofiber type composition analysis
Timing: Half a day (depends on the dataset size) In this section, we explain the filtering of non-myofiber ROIs and myofiber type classification. Here, we perform all the analyses in RStudio Software (v1.3.959 used for this protocol) 9 using R Statistical Software (v4.0.2 used for this protocol). 8 The R Markdown file is available on GitHub: https://github. com/tabbassidaloii/ImageProcessing/blob/main/MyofiberTyping/Rscript/MyofiberTyping.Rmd. In the steps below, we specify which R code chunk in this R markdown file should be used.
Note: The threshold used in the R markdown file is specific to one example dataset that was used in this protocol. The filtering threshold should be tuned for each dataset.

Exclusion of non-myofiber ROIs.
In image processing, the laminin segmentation is automated, which may include non-myofiber ROIs in the dataset. Therefore, the first step is to filter out the non-myofiber ROIs. We consider percentiles of 'segmentation certainty', CSA, and circularity values.
Note: Filtering out using a statistical value could result in removing true myofibers ROIs. However, due to the large myofiber ROIs, the dataset outcome is unaffected by small changes in percentile cutoff values.
Note: We verify the thresholds applied by visualizing the filtered ROIs using ''6.Visual_ Check_Filtering.ijm'' macro available on GitHub: https://github.com/tabbassidaloii/Image Processing/tree/main/MyofiberTyping/Macros. The visualization can be performed after each filtering step to justify the threshold. Here we only show an example of filtered myofibers after the last filtering step.
a. Filtering based on segmentation certainty. The measurement in the 'classification' channel can be used to assess the segmentation certainty and filter out non-myofiber ROIs ( Table 2).
Note: The measurements which are based on the ''Mean gray value'' are: Mean, laminin intensity inside the object (the smaller the better); Mean_boundary, laminin intensity on the boundary (the larger the better); StdDev_boundary, standard deviation of laminin intensity on the boundary (the smaller the better).
i. Pool all the data from all the samples by running ''readDataset'' R code chunk.
ii. Draw a density plot for each metric to define a filtering threshold by running ''Filt1_seg-mentationMetrics_denPlot'' R code chunk. iii. Define the threshold based on the density distribution and apply filtering for Mean and Mean_boundary by running ''Filt1_segmentationMetrics_filtering'' R chunk code. In our dataset, including 369073 ROIs, we excluded ROIs from the top and bottom 5th percentile (Figure 9).
Note: In this example, filtering based on StdDev_boundary is not used, because the StdDev_boundary has been improved by filtering for two other metrics.
Note: The filtering thresholds should be tuned for other input datasets. We recommend iteratively assessing the thresholds applied by visualizing the filtered ROIs. The Mean values for channels 1, 2 and 3 will be used to determine myofiber type in subsequent analysis. Channel 4 is not used in subsequent analysis. Channel 5 is the 'pixel classification score' that enables us to verify the quality of the segmentation. Please note that some values are duplicates (e.g., area and circularity only depend on the shape of the ROIs, and are therefore the same for all channels), and not all values are used for downstream analysis, these have been grayed out. Label consists of a unique identifier of the ROI and the channel number. Area is the cross-sectional area (CSA). Mean is the mean gray value intensity inside the object. Circ. is the circularity. Mean_boundary and StdDev_boundary are the mean gray value intensity and the standard deviation of the intensity on the object's boundary. Mean_distance is the object's mean distance to the edge of the section.

OPEN ACCESS
b. Filtering based on CSA (mm2). The next filtering step uses the CSA: large CSA values could represent multiple myofibers that were not individually segmented due to the absence of laminin staining or small CSA values might be either small non-myofiber cells or a space between contiguous myofibers. i. Draw a density plot to define a filtering threshold by running ''Filt2_CSA_denPlot'' R code chunk ( Figure 10A). ii. Apply filtering for CSA by running ''Filt2_CSA_filtering'' R chunk code.
Note: Include all ROIs including those that are filtered in the segmentation certainty filtering step.
Note: Depending on the dataset characteristic, CSA could have a diverse distribution across sample groups. Therefore, when the distribution is not the same, the filtering based on CSA should be separately applied to each group. However, applying filtering separately may introduce a bias, therefore these decisions about the study design require careful consideration. This is an important consideration to preserve biological differences between sample groups. We show an example for CSA across different muscles, we notice that the mean CSA differed between muscles, therefore the CSA-based filtering was made per muscle group.
Note: When the normal CSA range for a given sample group is known, filtering for CSA can be implemented by removing ROIs out of range.
Note: Since there is a right-skewed distribution, one may use different thresholds for both sides of the distribution ( Figure 10B).
Note: Verify the thresholds can be made by visualization of the filtered ROIs using the ''6.Visual_Check_Filtering.ijm'' macro available on GitHub.
c. Filtering based on circularity.
The circularity is an ROI measure, which ranges between 0 (elongated shape) and 1 (circle). This measure is used to filter out longitudinally sectioned myofibers and elongated ROIs that may represent multiple myofibers (mis-segmented ROIs). i. Draw a density plot of the entire dataset to define a filtering threshold by running ''Filt3_circularity_denPlot'' R code chunk ( Figure 11A).

OPEN ACCESS
ii. Apply filtering for circularity by running ''Filt3_circularity_filtering'' R chunk code. Based on density distribution, we include ROIs with Circularity > 1st percentile ( Figure 11B). Here, we include ROIs with 10 th percentile < CSA < 99 th percentile. In this protocol, we used human samples collected from six leg muscles, gracilis (GR), semitendinosus (ST), rectus femoris (RF), vastus lateralis (VL), vastus medialis (VM), and gastrocnemius lateralis (GL) muscles. We also included biopsies from the middle (STD) and distal (STM) end of the semitendinosus muscle to investigate differences within one muscle. For this dataset, we tested different thresholds. Visualizing the myofibers excluded, we realized that if we exclude more than 1% on the right side, we would filter out 'real' myofibers (visual base). The small nonmyofiber ROIs are not entirely excluded, but if we exclude more than 10% from the left (small ROIs), we would also filter out 'real' myofibers. Note: Apply filtering based on the percentile of all ROIs including those that are filtered in the previous filtering steps. 39. Visualizing and justifying the thresholds. In this step, we perform a spatial assessment of the excluded and included myofibers by visualizing the filtered ROIs using the ''6.Visual_Check_Filtering.ijm'' macro available on GitHub. a. Save a tab-delimited text file (with ''_Filt'' extension) in the ROI subdirectory by running ''save_visualizationInput'' R chunk code.
In each filtering step, the script adds a column with 0 and 1 values to the tab-delimited text file specifying which ROIs are included (1) or excluded (0). The script also adds an extra column with values between 0 and 1 to show the aggregated filtering results. These values will be used to give a distinct color to ROI excluded in each filtering step: i. 0: ROIs filtered based on segmentation certainty. ii. 0.4: ROIs included after filtering based on laminin segmentation certainty.
iii. 0.7: ROIs included after filtering based on CSA. iv. 1: ROIs included after filtering based on circularity. b. Open Fiji. c. Run ''6.Visual_Check_Filtering.ijm'' macro to visualize the filtered ROIs. d. Provide directory with tiff images and ROI subdirectory. e. For each sample, a jpeg file that ends with ''_check4'' will be saved in the image directory (Examples are in Figure 12).
40. Selecting one replicate per sample.
In this step, one section per sample is selected for further analysis. We recommend selecting the section with the highest number of myofibers after all filtering steps considering the quality of the sections. We consider only samples with a minimum of one hundred myofibers for downstream analysis.

classification of myofiber types.
The protocol for data-driven myofiber classification is detailed in Raz et al. 13

OPEN ACCESS
Note: The optimum bandwidth (h) is arbitrary and should be selected based on the dataset characteristic. We suggest using h values ranging from 0.01 to 0.05 to find the optimum value. In general, a lower value assigns myofibers to many small clusters, whereas a higher value would result in grouping the myofibers in larger clusters. We expect to have three to eight clusters as potentially biologically relevant.
ii. Remove the cluster with a low proportion of myofibers (typically < $2%).
iii. Draw boxplots to visualize the myofibers clusters (an example is in Figure 13). iv. Generate a spatial visualization of the clustering assignment per myofiber (an example is in Figure 14).

Biological interpretation.
The proportion of myofibers in each cluster can be then calculated per sample and can be used to compare myofiber type composition.
Note: The expression of MyHC isoforms describes contraction capacity. Changes in myofiber type composition were reported in pathological and physiological conditions. A quantitative assessment of myofiber type composition can help in understanding disease progression, 13 and of a therapeutic assessment. 16 The biological relevance of myofiber type composition clusters is not fully exploited. However, it allows recognizing hybrids, which are implicated in aging. 13 Myofiber type clusters differ between muscle groups in both mouse and human, 1,17,18 although not fully understood, it could suggest a relevance for muscle physiology and function.

EXPECTED OUTCOMES
This protocol allows segmentation of the whole muscle section, quantification of the three main MyHC isoforms, and myofiber typing using a semi-automated pipeline. The main user interaction is to check and potentially correct automatically generated tissue masks to exclude tissue folds or other potential artifacts from the analysis. Even though the pipeline is highly automated, each individual step can be inspected to provide additional insight and to alleviate any potential issues. The main outcome of the image processing steps is a text file with the mean fluorescent intensity of these isoforms as well as other properties of each individual segmented myofiber including cross-sectional area and circularity. The file also includes the measurement for the 'classification' channel allowing assessment of the myofiber 'segmentation certainty'. In the myofiber type composition analysis steps, following the filtering of non-myofiber ROIs, the MFI of three MyHC is used to cluster myofibers and to calculate the myofiber type composition in each sample.

QUANTIFICATION AND STATISTICAL ANALYSIS
In this protocol we show how to use text files to calculate and compare the myofiber type compositions between different sample groups in R. But other parameters (e.g., cross-sectional area) can also be used for further comparison and statistical analysis. In addition, the text file can be imported into any other data analysis software such as Python or Excel.

LIMITATIONS
The laminin segmentation is performed in an automated manner because manual segmentation is nearly undoable for whole tissue sections. Some mis-segmentation is unavoidable, but the majority of mis-segmentations are removed using the filter procedures described above.
In addition, we only applied this protocol on human and mouse samples (see Figure S1), and one may need to optimize it for other species. Finally, our image processing and quantification steps may need to be adjusted for other imaging equipment and image analysis software.

TROUBLESHOOTING
The authors are available for troubleshooting and advice. Figure 13. Myofiber types clustering example Each dot represents a myofiber. In the dataset used for this protocol, we found three main clusters. Each myofiber cluster has a major MyHC isoform: MyHC2A is the major isoform in Cluster 1 (red), MyHC1 is the major isoform in cluster 2 (blue) and MyHC2X is the major isoform in cluster 3 (green). Boxplots show the MFI of the corresponding MyHC isoforms in three different clusters.

Problem 1
When applying a gaussian blur to reduce any possible artifacts due to this binary mask, some biological signals might be removed.

Potential solution
In this protocol, we use a gaussian blur of 4 pixels (sigma = 4), but if it removes too many details the sigma can be adjusted in the ''2.Masked_Lamin.ijm'' macro (line 46).

Problem 2
After Immunostaining, there might be more staining artifacts on the edge of the tissue.

Potential solution
We provide a measure showing the distance of each segmented myofiber to the edge of the tissue in the text output file. This measurement can be used to filter out myofibers with staining artifacts close to the tissue borders.

RESOURCE AVAILABILITY
Users can report any issues with running this protocol to https://github.com/tabbassidaloii/ ImageProcessing/issues.

Lead contact
Further information and requests for resources and reagents should be directed to and will be fulfilled by the lead contact Lenard M. Voortman (l.m.voortman@lumc.nl).

Materials availability
This study did not generate new unique reagents.
Data and code availability A small subset of the dataset generated in the current study has been deposited on figshare: https:// doi.org/10.6084/m9.figshare.21324072.v2. All the image processing macros and R scripts are provided on GitHub: https://github.com/tabbassidaloii/ImageProcessing/tree/main/MyofiberTyping (https://doi.org/10.5281/zenodo.7466595). Myofibers assigned to the small clusters or filtered in the filtering steps are shown in gray. Panel B was generated by visualization of the clustering result using the ImageJ macro ''6.Visual_Check_Filtering.ijm'', followed by a colorization step using image editing software. Scale bar 500 mm.